LOP Manuscript Figures

files

  [1] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_004-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_005-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_006-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_011-dosage_plot.tsv" 
  [5] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_012-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_015-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_016-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_025-dosage_plot.tsv" 
  [9] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_026-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_027-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_028-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_029-dosage_plot.tsv" 
 [13] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_030-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_034-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_037-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_041-dosage_plot.tsv" 
 [17] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_042-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_043-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_044-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_045-dosage_plot.tsv" 
 [21] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_046-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_048-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_058-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_060-dosage_plot.tsv" 
 [25] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_062-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_064-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_065-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_072-dosage_plot.tsv" 
 [29] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_074-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_076-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_077-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_079-dosage_plot.tsv" 
 [33] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_081-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_086-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_088-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_090-dosage_plot.tsv" 
 [37] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_097-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_108-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_1157-dosage_plot.tsv" "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_116-dosage_plot.tsv" 
 [41] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_1160-dosage_plot.tsv" "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_117-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_1180-dosage_plot.tsv" "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_119-dosage_plot.tsv" 
 [45] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_1190-dosage_plot.tsv" "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_1192-dosage_plot.tsv" "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_1196-dosage_plot.tsv" "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_122-dosage_plot.tsv" 
 [49] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_123-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_124-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_126-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_128-dosage_plot.tsv" 
 [53] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_131-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_171-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_172-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_198-dosage_plot.tsv" 
 [57] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_204-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_205-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_219-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_220-dosage_plot.tsv" 
 [61] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_222-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_224-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_227-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_234-dosage_plot.tsv" 
 [65] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_236-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_238-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_244-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_247-dosage_plot.tsv" 
 [69] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_248-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_251-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_254-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_255-dosage_plot.tsv" 
 [73] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_256-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_257-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_258-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_259-dosage_plot.tsv" 
 [77] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_262-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_267-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_268-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_269-dosage_plot.tsv" 
 [81] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_270-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_272-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_273-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_274-dosage_plot.tsv" 
 [85] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_278-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_279-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_280-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_281-dosage_plot.tsv" 
 [89] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_286-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_289-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_292-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_293-dosage_plot.tsv" 
 [93] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_295-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_299-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_304-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_305-dosage_plot.tsv" 
 [97] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_306-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_308-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_309-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_311-dosage_plot.tsv" 
[101] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_313-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_349-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_361-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_363-dosage_plot.tsv" 
[105] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_365-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_372-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_373-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_387-dosage_plot.tsv" 
[109] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_396-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_397-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_398-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_401-dosage_plot.tsv" 
[113] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_403-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_406-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_407-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_414-dosage_plot.tsv" 
[117] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_415-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_417-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_420-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_421-dosage_plot.tsv" 
[121] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_425-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_428-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_429-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_430-dosage_plot.tsv" 
[125] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_433-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_435-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_439-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_449-dosage_plot.tsv" 
[129] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_450-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_452-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_454-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_456-dosage_plot.tsv" 
[133] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_459-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_460-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_461-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_463-dosage_plot.tsv" 
[137] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_464-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_465-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_466-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_467-dosage_plot.tsv" 
[141] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_469-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_473-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_474-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_477-dosage_plot.tsv" 
[145] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_479-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_480-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_483-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_485-dosage_plot.tsv" 
[149] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_487-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_488-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_493-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_494-dosage_plot.tsv" 
[153] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_497-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_499-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_500-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_502-dosage_plot.tsv" 
[157] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_505-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_509-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_510-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_513-dosage_plot.tsv" 
[161] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_514-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_522-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_523-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_525-dosage_plot.tsv" 
[165] "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_529-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_532-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_535-dosage_plot.tsv"  "../experiments/2_low_cov_cnv/data/plots/2x_LOP868_538-dosage_plot.tsv" 
[169] "../experiments/2_low_cov_cnv/data/plots/4x_LOP868-dosage_plot.tsv"     

names(all_dosage) 

  [1] "2x_LOP868_004"  "2x_LOP868_005"  "2x_LOP868_006"  "2x_LOP868_011"  "2x_LOP868_012"  "2x_LOP868_015"  "2x_LOP868_016"  "2x_LOP868_025"  "2x_LOP868_026"  "2x_LOP868_027"  "2x_LOP868_028"  "2x_LOP868_029"  "2x_LOP868_030"  "2x_LOP868_034"  "2x_LOP868_037"  "2x_LOP868_041" 
 [17] "2x_LOP868_042"  "2x_LOP868_043"  "2x_LOP868_044"  "2x_LOP868_045"  "2x_LOP868_046"  "2x_LOP868_048"  "2x_LOP868_058"  "2x_LOP868_060"  "2x_LOP868_062"  "2x_LOP868_064"  "2x_LOP868_065"  "2x_LOP868_072"  "2x_LOP868_074"  "2x_LOP868_076"  "2x_LOP868_077"  "2x_LOP868_079" 
 [33] "2x_LOP868_081"  "2x_LOP868_086"  "2x_LOP868_088"  "2x_LOP868_090"  "2x_LOP868_097"  "2x_LOP868_108"  "2x_LOP868_1157" "2x_LOP868_116"  "2x_LOP868_1160" "2x_LOP868_117"  "2x_LOP868_1180" "2x_LOP868_119"  "2x_LOP868_1190" "2x_LOP868_1192" "2x_LOP868_1196" "2x_LOP868_122" 
 [49] "2x_LOP868_123"  "2x_LOP868_124"  "2x_LOP868_126"  "2x_LOP868_128"  "2x_LOP868_131"  "2x_LOP868_171"  "2x_LOP868_172"  "2x_LOP868_198"  "2x_LOP868_204"  "2x_LOP868_205"  "2x_LOP868_219"  "2x_LOP868_220"  "2x_LOP868_222"  "2x_LOP868_224"  "2x_LOP868_227"  "2x_LOP868_234" 
 [65] "2x_LOP868_236"  "2x_LOP868_238"  "2x_LOP868_244"  "2x_LOP868_247"  "2x_LOP868_248"  "2x_LOP868_251"  "2x_LOP868_254"  "2x_LOP868_255"  "2x_LOP868_256"  "2x_LOP868_257"  "2x_LOP868_258"  "2x_LOP868_259"  "2x_LOP868_262"  "2x_LOP868_267"  "2x_LOP868_268"  "2x_LOP868_269" 
 [81] "2x_LOP868_270"  "2x_LOP868_272"  "2x_LOP868_273"  "2x_LOP868_274"  "2x_LOP868_278"  "2x_LOP868_279"  "2x_LOP868_280"  "2x_LOP868_281"  "2x_LOP868_286"  "2x_LOP868_289"  "2x_LOP868_292"  "2x_LOP868_293"  "2x_LOP868_295"  "2x_LOP868_299"  "2x_LOP868_304"  "2x_LOP868_305" 
 [97] "2x_LOP868_306"  "2x_LOP868_308"  "2x_LOP868_309"  "2x_LOP868_311"  "2x_LOP868_313"  "2x_LOP868_349"  "2x_LOP868_361"  "2x_LOP868_363"  "2x_LOP868_365"  "2x_LOP868_372"  "2x_LOP868_373"  "2x_LOP868_387"  "2x_LOP868_396"  "2x_LOP868_397"  "2x_LOP868_398"  "2x_LOP868_401" 
[113] "2x_LOP868_403"  "2x_LOP868_406"  "2x_LOP868_407"  "2x_LOP868_414"  "2x_LOP868_415"  "2x_LOP868_417"  "2x_LOP868_420"  "2x_LOP868_421"  "2x_LOP868_425"  "2x_LOP868_428"  "2x_LOP868_429"  "2x_LOP868_430"  "2x_LOP868_433"  "2x_LOP868_435"  "2x_LOP868_439"  "2x_LOP868_449" 
[129] "2x_LOP868_450"  "2x_LOP868_452"  "2x_LOP868_454"  "2x_LOP868_456"  "2x_LOP868_459"  "2x_LOP868_460"  "2x_LOP868_461"  "2x_LOP868_463"  "2x_LOP868_464"  "2x_LOP868_465"  "2x_LOP868_466"  "2x_LOP868_467"  "2x_LOP868_469"  "2x_LOP868_473"  "2x_LOP868_474"  "2x_LOP868_477" 
[145] "2x_LOP868_479"  "2x_LOP868_480"  "2x_LOP868_483"  "2x_LOP868_485"  "2x_LOP868_487"  "2x_LOP868_488"  "2x_LOP868_493"  "2x_LOP868_494"  "2x_LOP868_497"  "2x_LOP868_499"  "2x_LOP868_500"  "2x_LOP868_502"  "2x_LOP868_505"  "2x_LOP868_509"  "2x_LOP868_510"  "2x_LOP868_513" 
[161] "2x_LOP868_514"  "2x_LOP868_522"  "2x_LOP868_523"  "2x_LOP868_525"  "2x_LOP868_529"  "2x_LOP868_532"  "2x_LOP868_535"  "2x_LOP868_538"  "4x_LOP868"     

all_dosage <- files %>% 
  map(read_delim, delim = " ", col_names = T)
  
for (i in 1:length(files)) {
  all_dosage[[i]]$sample <- str_replace(files[i], "../experiments/2_low_cov_cnv/data/plots/", "")
  all_dosage[[i]]$sample <- gsub("-dosage_plot.tsv", "", all_dosage[[i]]$sample)
  all_dosage[[i]]$sample_num <- i
}

bin_dosage <- all_dosage %>% 
  bind_rows() %>% 
  mutate(sample = str_replace(sample, "-dosage_plot.tsv", "")) %>% 
  filter(sample != "2x_LOP868_279")

chrom_dosage <- bin_dosage %>% 
  group_by(chrom, sample, sample_num) %>% 
  summarize(chrom_readcount = sum(readcount)) %>% 
  arrange(sample, chrom) %>% 
  filter(!is.na(chrom)) %>% 
  ungroup()

with_totals <- chrom_dosage %>% 
  group_by(sample, sample_num) %>% 
  summarize(total_readcount = sum(chrom_readcount)) %>% 
  left_join(chrom_dosage, .) %>% 
  ungroup()

Joining, by = c("sample", "sample_num")

with_control <- with_totals %>% 
  filter(sample == "4x_LOP868") %>% 
  mutate(chromshort = str_replace(chrom, "chr", "")) %>% 
  mutate(chromshort = str_replace(chromshort, "^0", "")) %>% 
  rename(control_chrom_readcount = chrom_readcount) %>% 
  rename(control_total_readcount = total_readcount) %>%
  select(-sample, -sample_num) %>% 
  left_join(with_totals, .) %>% 
  mutate(normcov = 2* (chrom_readcount / total_readcount) / (control_chrom_readcount / control_total_readcount)) %>% 
  filter(!sample %in% c("4x_LOP868", "2x_LOP868_279"))

Joining, by = "chrom"

fig2A <- with_control %>% 
  filter(normcov <= mean(with_control$normcov + 3*sd(with_control$normcov))) %>% 
  ggplot(., aes(x = sample_num, y = normcov)) +
  geom_vline(xintercept = c(with_control$sample_num[which(with_control$normcov >= mean(with_control$normcov) + 3 * sd(with_control$normcov))]), color = "gray70") +
  geom_point() +
  geom_text(data = filter(with_control, normcov >= mean(with_control$normcov) + 3*sd(with_control$normcov)), aes(x = sample_num, y = normcov, label = chromshort), size = 5, fontface = "bold") +
  geom_line(y = mean(with_control$normcov), color = "green") +
  geom_line(y = mean(with_control$normcov) + 3 * sd(with_control$normcov), color = "red") +
  scale_x_continuous(breaks = c(1,167), labels = c("1", "167")) +
  scale_y_continuous(limits = c(1.6, 3.3)) +
  ggtitle("A") +
  labs(x = "Individual", y = "Chromosome\nCopy Number") +
  theme_bw() +
  theme(panel.grid = element_blank(),
        panel.background=element_rect(fill="white",color="black"),
        axis.title.x=element_text(size=18),
        axis.title.y=element_text(size=18,angle= 90, vjust=0.5), plot.title=element_text(size=24,face="bold",hjust=0,color="black"),
        axis.text.x=element_text(size=18,color="black"),
        axis.text.y=element_text(size=18,color="black"),
        panel.grid.major.x = element_line(color = "gray70"))
fig2A

# ggsave("Fig2A.pdf", plot=fig2A, width = 18, height = 3, units = "in", device = "pdf")

This is a mathematical representation of trisomy, so I can use this to go back and make overlay plots wihthout having to manually specify them.

trisomics <- with_control$sample[which(with_control$normcov >= mean(with_control$normcov) + 3 * sd(with_control$normcov))]

chrom_overlay <- function(df, chr) {
  plt <- df %>% 
    filter(chrom == chr) %>% 
    filter(sample %in% trisomics) %>% 
    filter(binsize >= 2.5e5) %>% 
    ggplot(., aes(x = start, y = normcov, fill = sample)) +
    geom_line(aes(color = sample), size = 1) +
    geom_line(data = filter(df, chrom == chr & !sample %in% trisomics & binsize >= 2.5e5), aes(x = start, y = normcov, fill = sample), alpha = 0.2, size = 0.5) +
    facet_wrap(~chrom, strip.position = "r") +
    scale_x_continuous(breaks = seq(0,8e7,by=2e7), labels = seq(0,80,by=20)) +
    guides(fill = F, color = F) +
    labs(x = "Position (Mb)", y = "Standardized\nCoverage") +
    theme_bw() +
      theme(
      panel.grid.minor.x = element_blank(),
      panel.grid.minor.y = element_blank(),
      panel.grid.major.y = element_line(color="gray30", size = 0.2),
      panel.grid.major.x = element_blank(),
      strip.background = element_rect(fill = "white", color = "black"),
      panel.background = element_rect(fill = "white", color = "black"),
      axis.text = element_text(size = 14),
      axis.title = element_text(size = 16),
      strip.text = element_text(size = 14))
  # return(plt)
  print(plt)
  # ggsave(paste(chr, "_dosage_overlay.pdf", plot = plt, width = 6, height = 2, units = "in", device = "pdf")) # save these if needed
}

fig3a <- chrom_overlay(bin_dosage, "chr01") + ggtitle("A")

Last part, do high coverage CNV here to finish Fig. S3

sc

[1]  80.67005  89.43656  71.43014 110.69950

fig_4d <- chrom_overlay(bin_dosage, "chr06") + ggtitle("D")

Ignoring unknown aesthetics: fill

fig_4d

fig4 <- grid.arrange(fig_4a, fig_4b, fig_4c, fig_4d, nrow = 4)
ggsave("Fig4.png", plot = fig4, width = 12, height = 12, units = "in", device = "png")
ggsave("Fig4.pdf", plot = fig4, width = 12, height = 12, units = "in", device = "pdf")

for (i in sprintf('chr%02d', 1:12)) {
  plt <- chrom_overlay(bin_dosage, i)
  dplt <- dp %>% 
  filter(chrom == i) %>% 
  filter(N_content <= 0.3) %>% 
  filter(sample == "4x_LOP868") %>% 
  ggplot(., aes(x = start, y = gc_norm_dp)) + 
  geom_point(alpha = 0.3, size = 0.8) +
  labs(x = "", y = "Alca Tarma\nCopy Number") +
  # ggtitle("A") +
  scale_y_continuous(limits = c(0, 220), breaks = c(0, sc[4]*.25, sc[4]*.5, sc[4]*0.75, sc[4], sc[4]*1.25, sc[4]*1.5, sc[4]*1.75, sc[4]*2), labels = c(0,1,2,3,4,5,6,7,8)) +
  scale_x_continuous(limits = c(0,6e7), breaks = seq(0, 6e7, by = 2e7), labels = seq(0, 60, by = 20)) +
  facet_wrap(~chrom, nrow = 12, strip.position = "right") +
  theme(
    panel.grid.major.x = element_blank(),
    panel.grid.minor.x = element_blank(),
    panel.grid.major.y = element_line(color = "gray30", size = 0.2),
    panel.grid.minor.y = element_blank(),
    panel.background = element_rect(fill = "white", color = "black"),
    axis.text = element_text(size = 14),
    axis.title = element_text(size = 16),
    strip.text = element_text(size = 14),
    # plot.title = element_text(size = 24, face = "bold"),
    strip.background = element_rect(fill = "white", color = "black")
  )
  comb <- grid.arrange(dplt,plt)
  ggsave(paste0(i,"_figS3_combo_plot.png"), plot = comb, width = 9, height = 4, units = "in", device = "png")
}

Ignoring unknown aesthetics: fill

Fig. S5

names(sc) <- str_extract(dp_files, "LOP.+") %>%
  str_replace(., ".tsv", "") %>% 
  str_replace(., "_", ".")
names(sc)

overlay_emphasize_one <- function(df, dihap, chr) {
  plt <- df %>%
    filter(chrom == chr) %>% 
    filter(sample != dihap) %>% 
    filter(binsize >= 2.5e5) %>%
    ggplot(., aes(x = start, y = normcov, fill = sample)) +
    geom_line(alpha = 0.2) +
    geom_line(data = filter(df, chrom == chr & sample == dihap & binsize >= 2.5e5), color = "red", size = 1.2) +
    scale_x_continuous(breaks = seq(0,6e7,by=2e7), labels=seq(0,60,by=20)) +
    facet_wrap(~chrom, strip.position = "r") +
    labs(x = "Position (Mb)", y = "Standardized Coverage") +
    theme_bw() +
    theme(strip.background = element_rect(fill = "white", color = "black"),
          axis.text = element_text(size = 14),
          axis.title = element_text(size = 16),
          strip.text = element_text(size = 14),
          panel.grid.major.x = element_blank(),
          panel.grid.minor = element_blank())
  return(plt)
}
hicov <- function(df, dihap, chr, s, e) {
  
  sc_to_use <- dihap %>% 
    str_remove("[24]x_") %>% 
    str_replace(., "_", ".")
  # print(sc_to_use)
  
  sc_loc <- sc[which(names(sc) == sc_to_use)]
  
  plt <- df %>% 
    filter(chrom == chr) %>% 
    filter(between(start, s, e)) %>% 
    filter(N_content <= 0.3) %>% 
    filter(sample == dihap) %>% 
    ggplot(., aes(x = start, y = gc_norm_dp)) +
    # geom_point() +
    # geom_line() +
    labs(x = "Position (Mb)", y = paste(sc_to_use, "Copy Number", sep = " ")) +
    facet_wrap(~chrom, strip.position = "r") +
    scale_x_continuous(breaks = seq(s,e, by = 2.5e5), labels = seq(s/1e6, e/1e6, by = 0.25)) +
    scale_y_continuous(limits = c(0, 2*sc_loc), breaks = c(0, 0.5*sc_loc, sc_loc, 1.5*sc_loc, 2*sc_loc), labels = seq(0,4)) +
    theme_bw() +
    theme(strip.background = element_rect(fill = "white", color = "black"),
          axis.text = element_text(size = 14),
          axis.title = element_text(size = 16),
          strip.text = element_text(size = 14),
          panel.grid.major.x = element_blank(),
          panel.grid.minor = element_blank())
  return(plt)
}

s5_1 <- overlay_emphasize_one(bin_dosage, "2x_LOP868_004", "chr06")
s5_2 <- overlay_emphasize_one(bin_dosage, "2x_LOP868_064", "chr06")
s5_3 <- overlay_emphasize_one(bin_dosage, "2x_LOP868_305", "chr06")
s5_4 <- hicov(dp, "2x_LOP868_004", "chr06", 3.2e7, 3.4e7) + geom_point() + geom_line()
s5_5 <- hicov(dp, "2x_LOP868_064", "chr06", 3.2e7, 3.4e7) + geom_point() + geom_line()
s5_6 <- hicov(dp, "2x_LOP868_305", "chr06", 3.2e7, 3.4e7) + geom_point() + geom_line()
s5_7 <- hicov(dp, "2x_LOP868_004", "chr06", 0, 6e7) + geom_point(alpha = 0.3, size = 0.8) + scale_x_continuous(limits = c(0, 6e7), breaks = seq(0, 6e7, by = 2e7), labels = seq(0, 60, by = 20))

Scale for 'x' is already present. Adding another scale for 'x', which will replace the existing scale.

s5_8 <- hicov(dp, "2x_LOP868_064", "chr06", 0, 6e7) + geom_point(alpha = 0.3, size = 0.8) + scale_x_continuous(limits = c(0, 6e7), breaks = seq(0, 6e7, by = 2e7), labels = seq(0, 60, by = 20))

Scale for 'x' is already present. Adding another scale for 'x', which will replace the existing scale.

s5_9 <- hicov(dp, "2x_LOP868_305", "chr06", 0, 6e7) + geom_point(alpha = 0.3, size = 0.8) + scale_x_continuous(limits = c(0, 6e7), breaks = seq(0, 6e7, by = 2e7), labels = seq(0, 60, by = 20))

Scale for 'x' is already present. Adding another scale for 'x', which will replace the existing scale.

Fig S4: Kernel density estimator plots illlustrating dosage variant clustering to get genotypes Note, red points were drawn in as manual points, using x values specified according to the cluster center, and the Y axis point drawn s/t it was at the peak of the kernel density function.

kde_bin <- function(df, chr, s) {
  active_bin <- filter(df, chrom == chr & start == s) %>% 
    mutate(bin = paste(chrom, start, end, sep = "_")) 
  plt <- ggplot(active_bin, aes(x = normcov)) +
    geom_histogram(aes(y = ..density..), binwidth = 0.1, fill = "blue", alpha = 0.5) +
    geom_density(kernel = "ep", bw = quantile(dist(active_bin$normcov), 0.4)) +
    scale_x_continuous(limits = c(0,4)) +
    facet_wrap(~bin, strip.position = "r")
  return(plt)
}

s4_1 <- kde_bin(bin_dosage, "chr01", 4.25e6)
s4_2 <- kde_bin(bin_dosage, "chr01", 2.425e7)
s4_3 <- kde_bin(bin_dosage, "chr01", 7.95e7)
ggsave("Fig_S4_1.png", plot = s4_1, width = 4, height = 2, units = "in", device = "png")
ggsave("Fig_S4_2.png", plot = s4_2, width = 4, height = 2, units = "in", device = "png")
ggsave("Fig_S4_3.png", plot = s4_3, width = 4, height = 2, units = "in", device = "png")

Fig. S6: Power analysis of rare indels in LOP dihaploid population. Note, I added the figure legend in Affinity Designer by copying the legend from figure s6_legend

simhyb_les_alleles <- read_tsv("../experiments/3_fake_hybrids/indel_bin_alleles_simhyb.txt", col_names = T) %>% 
  mutate(control_group = ifelse(grepl("SRR6123031", sample), "Simulated Alca Tarma x PL4",
                                ifelse(grepl("SRR6123183", sample), "Simulated Alca Tarma x IVP101", "Simulated Matrilineal"))) %>% 
  mutate(control_group = factor(control_group, levels = c("Simulated Matrilineal", "Simulated Alca Tarma x IVP101", "Simulated Alca Tarma x PL4"))) %>% 
  mutate(bin = paste(chrom, as.integer(start), as.integer(end), sep = "_"))

Parsed with column specification:
cols(
  chrom = col_character(),
  start = col_double(),
  end = col_double(),
  size = col_double(),
  sample = col_character(),
  PercALesion = col_double(),
  TotalCovLesion = col_double(),
  NumSNPLesion = col_double(),
  PercANotLesion = col_double(),
  TotalCovNotLesion = col_double(),
  NumSNPNotLesion = col_double()
)

head(simhyb_les_alleles)

plot_power_simhyb_lesion <- function(df, chr, s) {
  plt <- df %>% 
    filter(chrom == chr) %>% 
    filter(start == s) %>% 
    ggplot(., aes(x = PercALesion, fill = control_group)) +
    geom_histogram(binwidth = 0.01, alpha = 0.8) +
    labs(x = "Percent HI Allele", y = "Count") +
    facet_wrap(~bin) +
    theme_bw() +
    # guides(fill = F) +
    scale_x_continuous(limits = c(-0.01,1)) +
    scale_fill_manual(values = c("#00BA38", "#F8766D", "#00BFC4")) +
    theme(panel.grid = element_blank(),
          axis.text = element_text(size = 14),
          axis.title = element_text(size = 16),
          strip.text = element_text(size = 14),
          strip.background = element_rect(fill = "white", color = "black"))
  return(plt)
}

Fig. S2: Binned SNP result, including simulated hybrid power test.

melt_bin_alleles <- function(x) {
  df <- read_tsv(x, na = ".")
  
  ObsHI <- df %>% 
    select(Chrom,Start,End,Max,matches("Obs%A$")) %>% 
    mutate(chrbin = floor(Start/df$End[1])) %>% 
    gather(Ind, ObsPerHI, -Chrom, -Start, -End, -Max, -chrbin) %>% 
    mutate(Ind = str_replace_all(Ind, "v-(.+)-Obs%A", "\\1"))
  
  CalcHI <- df %>% 
    select(Chrom, Start, End, Max, matches('Calc%A$')) %>% 
    mutate(chrbin = floor(Start/df$End[1])) %>% 
    gather(Ind, CalcPerHI, -Chrom, -Start, -End, -Max, -chrbin) %>% 
    mutate(Ind = str_replace_all(Ind, "v-(.+)-Calc%A", "\\1")) 
  
  Cov <- df %>% 
    select(Chrom, Start, End, Max, matches('Cov$')) %>% 
    mutate(chrbin = floor(Start/df$End[1])) %>% 
    gather(Ind, Cov, -Chrom, -Start, -End, -Max, -chrbin) %>% 
    mutate(Ind = str_replace_all(Ind, "v-(.+)-Cov", "\\1"))
  
  df2 <- ObsHI %>% 
    inner_join(.,Cov) %>% 
    mutate(HIcalls = round(Cov * ObsPerHI / 100)) %>% 
    filter(Chrom != "chr00") %>% 
    # filter(!grepl("sub", Ind)) %>% 
    mutate(prenducer = str_extract(Ind, "(SRR6123031)|(SRR6123183)")) %>% 
    mutate(Inducer = ifelse(prenducer == "SRR6123031", "PL4", NA)) %>% 
    mutate(Inducer = ifelse(prenducer == "SRR6123183", "IVP101", Inducer)) %>% 
    mutate(ObsPerHI_filt = ifelse(Cov >= 10, ObsPerHI, NA)) %>% 
    mutate(start_cen = Start + (End-Start)/2) 
  print(head(df2))
  
  return(df2)
}

big <- melt_bin_alleles("../experiments/3_fake_hybrids/1Mb_bin_alleles_simhyb.txt") %>% 
  mutate(`Control Group` = ifelse(Inducer == "PL4", "Simulated Alca Tarma x PL4", ifelse(Inducer == "IVP101", "Simulated Alca Tarma x IVP101", "Simulated Matrilineal")))

Parsed with column specification:
cols(
  .default = col_double(),
  Chrom = col_character()
)
See spec(...) for full column specifications.
Joining, by = c("Chrom", "Start", "End", "Max", "chrbin", "Ind")

pc <- big %>% 
  filter(grepl("SRR6123031|SRR6123183", Ind))
nc <- big %>% 
  filter(grepl("uniparental", Ind))

# Define bins
bigp <- big %>% 
  select(Chrom, Start, End) %>% 
  distinct() %>% 
  mutate(overlap_test = NA) %>% 
  mutate(missing_n_test = NA) %>%
  mutate(missing_n_frac = NA) %>% 
  mutate(missing_p_test = NA) %>% 
  mutate(missing_p_frac = NA)
nrow(bigp)

[1] 730

# define critical fraction of missing data
nafrac_crit <- 0.05 # may have to futz with this

for (i in 1:nrow(bigp)) {
  
  chrom <- bigp$Chrom[i]
  strt <- bigp$Start[i]
  endd <- bigp$End[i]
  
  # Define range of values for bin in negative control population
  rng_p <- pc %>% 
    filter(Chrom == chrom) %>% 
    filter(Start == strt) %>% 
    filter(End == endd) %>% 
    select(ObsPerHI)
  
  pmin <- min(rng_p$ObsPerHI, na.rm = T)
  print(paste("pmin:", pmin))
  
  pmax <- max(rng_p$ObsPerHI, na.rm = T)
  print(paste("pmax:", pmax))
  
  nafrac_p <- sum(is.na(rng_p$ObsPerHI)) / length(rng_p$ObsPerHI)
  
  # Define range of values for bin in positive control population
  rng_n <- nc %>% 
    filter(Chrom == chrom) %>% 
    filter(Start == strt) %>% 
    filter(End == endd) %>% 
    select(ObsPerHI)
  
  nmin <- min(rng_n$ObsPerHI, na.rm = T)
  print(paste("nmin:", nmin))
  
  nmax <- max(rng_n$ObsPerHI, na.rm = T)
  print(paste("nmax:", nmax))
  
  nafrac_n <- sum(is.na(rng_n$ObsPerHI)) / length(rng_n$ObsPerHI)
  
  # Define PASS/FAIL criteria as variables in original data frame
  # first failure criterion: if max of negative control greater than min of positive control
  bigp$overlap_test[i] = ifelse(nmax >= pmin, "FAIL", "PASS") # can then cbind when I'm done.
  
  # second failture criterion: missing data in negative control
  bigp$missing_n_test[i] <- ifelse(nafrac_n >= nafrac_crit, "FAIL", "PASS")
  bigp$missing_n_frac[i] <- nafrac_n
  
  # third failure criterion: missing data in positive control
  bigp$missing_p_test[i] <- ifelse(nafrac_p >= nafrac_crit, "FAIL", "PASS")
  bigp$missing_p_frac[i] <- nafrac_p
}

qc_bins <- left_join(pc, bigp) %>% 
  mutate(ObsPerHI_qc = ifelse(overlap_test == "PASS", ObsPerHI, NA))

Joining, by = c("Chrom", "Start", "End")

real <- melt_bin_alleles("../experiments/4_low_cov_snp/1Mb_bin_alleles_LOP.txt")

Parsed with column specification:
cols(
  .default = col_double(),
  Chrom = col_character()
)
See spec(...) for full column specifications.
Joining, by = c("Chrom", "Start", "End", "Max", "chrbin", "Ind")

real3 <- real %>% 
  filter(!(Ind %in% c("2x_IVP_101", "2x_PL_4", "4x_LOP868", "2x_LOP868_279"))) %>% 
  filter(Chrom != "chr00") %>% 
  left_join(., bigp) %>% 
  mutate(ObsPerHI_qc = ifelse(overlap_test == "PASS", ObsPerHI, NA)) %>% 
  mutate(ObsPerHI_qc_cov = ifelse(Cov < 30, NA, ObsPerHI_qc)) %>% 
  mutate(ObsPerHI_qc_50cov = ifelse(Cov < 50, NA, ObsPerHI_qc))

Joining, by = c("Chrom", "Start", "End")

fig_s2 <- real3 %>% 
  ggplot(., aes(x = Start, y = ObsPerHI_qc_cov, fill = Ind)) +
  geom_line(size = 0.5, alpha = 0.5) +
  geom_line(data = qc_bins, aes(x = Start, y = ObsPerHI_qc, fill = Ind, color = `Control Group`), size = 0.2, alpha = 0.3) +
  facet_wrap(~Chrom, nrow = 6, strip.position = "r") +
  labs(x = "Position (Mb)", y = "%HI SNPs") +
  scale_x_continuous(limits = c(0,8e7), breaks = seq(0,8e7,by=2e7), labels=seq(0,80,by=20)) +
  scale_y_continuous(limits = c(0,100)) +
  # scale_color_manual(values = c("#00BA38", "#619CFF", "#F8766D")) +
  theme_bw() +
  theme(panel.grid.minor = element_blank(),
        panel.grid.major.x = element_blank(),
        panel.grid.major.y = element_line(linetype = "dotted"),
        legend.position = "bottom",
        axis.text = element_text(size = 14),
        axis.title = element_text(size = 16),
        strip.text = element_text(size = 14))

Ignoring unknown aesthetics: fill

fig_s2
ggsave("Fig_S2.png", width = 10, height = 10, units = "in", device = "png", plot = fig_s2)

ggsave("Fig_S2.pdf", width = 10, height = 10, units = "in", device = "pdf", plot = fig_s2)

Overlay combined dosage and 1Mb SNP plots for Fig. S1. These were later put together in Affinity Designer with metaphase spreads, if spreads were available. TODO: put raw images of spreads in a separate folder.

numblanks_snp <- 10
mid <- function(x) {
  n = floor(nrow(x)/2)
  return(x$bin[n])
}

plot.snp <- function(x,y,z) {
  keep <- real3 %>% 
    filter(Ind == x) %>% 
    select(Chrom, Start, End, Max, chrbin, Ind, ObsPerHI, Cov, HIcalls, ObsPerHI_qc_cov)
    
  keep$bin <- seq(1:nrow(keep))
  # print(head(keep))
  
  stuf <- c(rep(NA, numblanks_snp))
  stuffer <- data.frame("Chrom"=stuf,
                        "Start"=stuf,
                        "End"=stuf,
                        "Max"=stuf,
                        "chrbin"=stuf,
                        "Ind"=stuf,
                        "ObsPerHI"=stuf,
                        "Cov"=stuf,
                        "HIcalls"=stuf,
                        "ObsPerHI_qc_cov"=stuf,
                        "bin"=stuf)
  chr.list <- split(keep, f=keep$Chrom)
  chr.list.stuffed <- lapply(chr.list[1:11], function(x) rbind(x,stuffer))
  ind.mod <- bind_rows(chr.list.stuffed, chr.list[[12]])
  ind.mod$bin2 <- seq(1:nrow(ind.mod))
  
  midpoints <- sapply(chr.list[1:12], mid)
  if (z == T) {
    x.axis.breaks <- which(ind.mod$bin %in% midpoints)
    x.axis.labels <- names(midpoints)
  } else {
    x.axis.breaks <- seq(1,12)
    x.axis.labels <- rep("",12)
  }
    plt.ind.JMPlike <- ggplot(ind.mod, aes(x=bin2,y=ObsPerHI_qc_cov)) +
    labs(x="", y="% HI SNPs") +
    ggtitle(" ") +
    # geom_line(y=0.0, color="black", linetype="dashed", size = 2) +
    # geom_line(y = 33/4, color = "black", linetype = "dashed") +
    geom_line(color=y, fill=NULL,size=3) +
    scale_color_manual(values=c(rep(y,12))) +
    theme(panel.grid.minor.x = element_blank(),
          panel.grid.minor.y=element_blank(),
          panel.grid.major.x = element_blank(),
          panel.grid.major.y = element_line(color="black", linetype="dashed"),
          panel.background = element_rect(fill="white",color="black"),
          axis.text.x = element_text(size=18, color="black"),
          axis.text.y=element_text(size=18,color="black"),
          axis.title.y=element_text(size=18, angle=90, vjust = 0.5, margin = margin(t = 0, r = 25, b = 0, l = 0)),
          axis.ticks = element_blank(),
          plot.title=element_text(size=24,face="bold",hjust=0),
          panel.border = element_rect(color = "black", fill = NA, size = 2)) +
    geom_point(size = 2.0,color="black") +
    guides(fill=FALSE, color=FALSE) +
    scale_y_continuous(limits=c(-5,100), breaks=c(100,66,33,0)) +
    scale_x_continuous(breaks=x.axis.breaks, labels=x.axis.labels)
  plt.ind.JMPlike
}

# retrieves midpoints of each chromosome
mid.dosage <- function(x) {
  n=floor(nrow(x)/2)
  return(x$bin[n])
}
# Dosage plot function arguments:
# x is the name of the .tsv fle to be read in
# y is the title of the plot. I gsub the string x to get the title. I could hard code this but I leave it as an arg to be more customizable.
plot.dosage <- function(ind, y) {
  todo <- which(names(all_dosage) == ind)
  # print(todo)
  
  plt.dosage <- all_dosage[[todo]] %>%
    ggplot(., aes(x = bin2, y = normcov)) +
    labs(x = "", y="Standardized\nCopy Number") +
    geom_line(color = "#008080", fill = "white", size = 3) +
    ggtitle(y) +
    geom_point(size = 2, color = "black") +
    guides(fill = F, color = F) +
    scale_y_continuous(limits = c(0,5), breaks = seq(0,5), labels = seq(0,5)) +
    theme(panel.grid.minor.x = element_blank(),
          panel.grid.minor.y=element_blank(),
          panel.grid.major.x = element_blank(),
          panel.grid.major.y = element_line(color="black", linetype="dashed"),
          panel.background = element_rect(fill="white",color="black"),
          axis.text.x = element_blank(),
          axis.text.y=element_text(size=18,color="black"),
          axis.title.y=element_text(size=18, angle=90, vjust = 0.5, margin = margin(t = 0, r = 25, b = 0, l = 0)),
          axis.ticks = element_blank(),
          plot.title=element_text(size=24,face="bold",hjust=0),
          panel.border = element_rect(fill = NA, color = "black", size = 2))
  plt.dosage
}

plot.dosage("2x_LOP868_064", "LOP868.064")

Ignoring unknown parameters: fill

plot.snp("2x_LOP868_004", "#008080", TRUE)

Ignoring unknown parameters: fill

stack.plot <- function(x) {
  # dosage.fh <- gsub("$", "-dosage_plot.tsv", x)
  # dosage <- read.table(dosage.fh, header=T)
  gt <- ggplotGrob(plot.dosage(x, ""))
  gb <- ggplotGrob(plot.snp(x, "#008080", TRUE))
  fig.s1.template <- rbind(gt,gb, size = "first")
  fig.s1.template$widths <- unit.pmax(gt$widths)
  fig.s1.template$layout[grepl("guide", fig.s1.template$layout$name), c("t", "b")] <- c(1, nrow(fig.s1.template))
  grid.newpage()
  grid.draw(fig.s1.template)
  # out.fh <- gsub("^", "2018_0531_", i)
  out.fh <- gsub("$", "_stack_plot.eps", x)
  ggsave(out.fh, plot=fig.s1.template, width = 24, height = 6, units = "in", device = "eps") # 18 inch width for supplemental figures, 24 inch width for main
}

for (i in names(all_dosage)) {
  stack.plot(i)
}

Ignoring unknown parameters: fillRemoved 497 rows containing missing values (geom_point).Ignoring unknown parameters: fillRemoved 1 rows containing missing values (geom_path).Removed 249 rows containing missing values (geom_point).

Circos plot generation done in a subdirectory of results